home *** CD-ROM | disk | FTP | other *** search
/ C/C++ Interactive Reference Guide / C-C++ Interactive Reference Guide.iso / c_ref / csource4 / 247_01 / bnprime.c < prev    next >
Text File  |  1989-04-19  |  4KB  |  173 lines
  1. /*
  2.  *   MIRACL prime number routines - test for and generate prime numbers
  3.  *   bnprime.c
  4.  */
  5.  
  6. #include <stdio.h>
  7. #include "miracl.h"
  8.  
  9. /* access global variables */
  10.  
  11. extern int depth;    /* error tracing ..*/
  12. extern int trace[];  /* .. mechanism    */
  13.  
  14. extern big w1;       /* workspace variables */
  15. extern big w2;
  16. extern big w3;
  17. extern big w4;
  18. extern big w5;
  19. extern big w6;
  20. extern big w7;
  21.  
  22.  
  23. bool prime(x)
  24. big x;
  25. {  /*  test for primality (probably) ;TRUE if x is  *
  26.     *  prime. test done NTRY times; chance of wrong * 
  27.     *  identification << (1/4)^NTRY                 */
  28.     int i,j,k,m,n,sx,nits;
  29.     static char primes[]={2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,
  30.                            61,67,71,73,79,83,89,97,101,103};
  31.     static long pmult[]={223092870L,58642669L,600662303L,33984931L,89809099L};
  32.     if (ERNUM) return TRUE;
  33.     sx=size(x);
  34.     if (sx<=1) return FALSE;
  35.     if (sx<105) 
  36.     {
  37.         for (i=0;i<27;i++)
  38.               if (sx==(int)primes[i]) return TRUE;
  39.         return FALSE;
  40.     }
  41.     depth++;
  42.     trace[depth]=22;
  43.     if (TRACER) track();
  44.     for (i=0;i<4;i++)
  45.     { /* check if divisible by first few primes */
  46.         lconv(pmult[i],w3); 
  47.         if (gcd(w3,x,w4)!=1)
  48.         { /* factor found... */
  49.             depth--;
  50.             return FALSE;
  51.         }
  52.     }
  53.     nits=logb2(x)/5;        /* lot less work than a powmod */
  54.     convert(5,w1);
  55.     convert(2,w2);
  56.     convert(1,w3);
  57.     k=m=1;
  58.     for (i=1;i<nits;i++)
  59.     { /* try a few iterations of Pollard-Rho factoring method */
  60.         k--;
  61.         if (k==0)
  62.         {
  63.             copy(w1,w2);
  64.             m*=2;
  65.             k=m;
  66.         }
  67.         convert(1,w4);
  68.         mad(w1,w1,w4,x,x,w1);
  69.         subtract(w2,w1,w4);
  70.         mad(w4,w3,w3,x,x,w3);
  71.     }
  72.     if (gcd(w3,x,w4)!=1)
  73.     { /* factor found... */
  74.         depth--;
  75.         return FALSE;
  76.     }
  77.     decr(x,1,w1); /* calculate k and w3 ...    */
  78.     k=0;
  79.     while (subdiv(w1,2,w1)==0)
  80.     {
  81.         k++;
  82.         copy(w1,w3);
  83.     }              /* ... such that x=1+w3*2**k */
  84.     for (n=0;n<NTRY;n++)
  85.     { /* perform test NTRY times */
  86.         convert((int)primes[n],w4);
  87.         j=0;
  88.         powmod(w4,w3,x,w4);
  89.         decr(x,1,w2);
  90.         while ((j>0 || size(w4)!=1) 
  91.               && compare(w4,w2)!=0)
  92.         {
  93.             j++;
  94.             if ((j>1 && size(w4)==1) || j==k)
  95.             { /* definitely not prime */
  96.                 depth--;
  97.                 return FALSE;
  98.             }
  99.             mad(w4,w4,w4,x,x,w4);
  100.         }
  101.     }
  102.     depth--;
  103.     return TRUE;  /* probably prime */
  104. }
  105.  
  106. void nxprime(w,x)
  107. big w,x;
  108. {  /*  find next highest prime from w using     * 
  109.     *  probabilistic primality test NTRY times  */
  110.     if (ERNUM) return;
  111.     depth++;
  112.     trace[depth]=21;
  113.     if (TRACER) track();
  114.     copy(w,x);
  115.     if (size(x)<2) 
  116.     {
  117.         convert(2,x);
  118.         depth--;
  119.         return;
  120.     }
  121.     if (subdiv(x,2,w1)==0) incr(x,1,x);
  122.     else                   incr(x,2,x);
  123.     while (!prime(x)) incr(x,2,x);
  124.     depth--;
  125. }
  126.  
  127. void strongp(p,n)
  128. big p;
  129. int n;
  130. { /* generate strong prime number p suitable *
  131.    * for encryption. p will be > n decimal   *
  132.    * digits in length. See Gordon J. 'Strong *
  133.    * primes are easy to find', Advances in   *
  134.    * Cryptology, Proc Eurocrypt 84, Lecture  *
  135.    * notes in Computer Science, Vol. 209     */
  136.     int k,n1,n2;
  137.     if (ERNUM) return;
  138.     depth++;
  139.     trace[depth]=47;
  140.     if (TRACER) track();
  141.     n2=n/3;
  142.     n1=2*n2+n%3;
  143.     bigdig(w6,n1,10);
  144.     nxprime(w6,w6);
  145.     bigdig(w5,n2,10);
  146.     nxprime(w5,w5);  /* find w5 and w6 prime */
  147.     k=0;
  148.     do
  149.     { /* find prime p = k.w6 + 1  with k even */
  150.         k+=2;
  151.         premult(w6,k,p);
  152.         incr(p,1,p);
  153.     } while (!ERNUM && !prime(p));
  154.     multiply(p,w5,w7);
  155.     decr(p,1,p);
  156.     powmod(w5,p,w7,w6);
  157.     incr(p,1,p);
  158.     decr(w5,1,w5);
  159.     powmod(p,w5,w7,w5);
  160.     subtract(w6,w5,w6);     /* w6 = (w5^(p-1) - p^(w5-1)) mod w7  */
  161.     if (subdiv(w6,2,w5)==0) add(w6,w7,w6);
  162.     k=0;
  163.     do
  164.     { /* find prime p = w6 + k.w7  with k even and p=2 mod 3 */
  165.         k+=2;
  166.         if (ERNUM) break;
  167.         premult(w7,k,p);
  168.         add(p,w6,p);
  169.     } while (subdiv(p,3,w5)!=2 || !prime(p));
  170.     depth--;
  171. }
  172.  
  173.